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We formulate a continuum linear response theory on the basis of the Hartree-Fock- 
Bogoliubov formalism in the coordinate space representation in order to describe low-lying 
and high-lying collective excitations which couple to one-particle and two-particle continuum 
states. Numerical analysis is done for the neutron drip-line nucleus 24 0. A low-lying collective 
mode that emerges above the continuum threshold with large neutron strength is analyzed. 
The collective state is sensitive to the density-dependence of the pairing. The present theory 
satisfies accurately the energy weighted sum rule. This is guaranteed by treating the pairing 
selfconsistently both in the static HFB and in the dynamical linear response equation. 

I. INTRODUCTION 

0\ 

Collective excitation in exotic unstable nuclei, especially in neutron-rich nuclei near the drip-line, is one of 
the most interesting nuclear structure issues. The presence of the neutron halo or skin structures, or more 
I/-) ' generally of loosely bound neutrons and the very shallow Fermi energy will modify characters of collective 
excitations known in stable nuclei. It is further interesting if new kinds of collective mode emerge. The 
linear response theory or the random phase approximation is one of the powerful tools to study such issues. 
Since the method itself is a general framework to describe normal modes of excitation built on a reference 
state given by mean-field approximations Q, it is advantageous to, rather than very light drip-line nuclei 
such as Li, heavier systems which will be studied in future experiments. Previous works in this direction 
have analyzed the giant resonances and the threshold neutron strength in exotic nuclei with closed shell 
configurations [^| [f] with use of the continuum linear response theory formulated in the coordinate space 
|aj|. However, to explore more systematically nuclei with open shell configurations, the theory has to be 
£j ' extended so as to take into account the pairing correlation, which may play essential roles especially for 
^ , low-lying collective excitations. 

Indeed the pairing correlation is a key element in the study of nuclei near the drip- line |7]-||]. A special 
feature of the pairing in drip-line nuclei is that the correlation takes place not only among bound orbits 
in the potential well but also in continuum states above the zero energy threshold. The Hartree-Fock- 
Bogoliubov (HFB) theory formulated in the coordinate space representation MM has a great power in this 
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respect since it allows one to treat properly the pairing in the continuum orbits, whereas the conventional 
BCS approximation is not suited for this purpose. The coordinate space HFB has been applied extensively 
and clarifying some aspects of pairing effects on the ground state properties including the halo and the 
skin |^|]|l(]-[l7j . Descriptions of deformed exotic nuclei are also under current developments 10-^] • The 
coordinate space HFB has been also used together with the relativistic mean-field models |23|-|25|| . 

The shallow Fermi energy in drip-line nuclei makes the threshold for particle continuum very low. It is 
therefore important to include effects of the continuum states in describing not only the ground state but also 
the excitation modes. Attempts have been made to describe both the pairing correlation and the continuum 
effects on collective states in the linear response formalism. These continuum quasiparticle linear response 
theories [^6|-[26[|, however, rely on the conventional BCS approximation, which may not be very accurate 
near the drip-lines. On the other hand, an approach is proposed to build a quasiparticle random phase 
approximation (QRPA) with use of the canonical basis in the coordinate space HFB j^o). However, the 
continuum effect is not precisely accounted in this approach since the single-quasiparticle basis is discretized. 
Other QRPA models using the discrete BCS quasiparticle basis f|l] |3j| have a similar problem. 

In the present paper, we extend the approach of Ref. and formulate a new continuum quasiparticle 
linear response theory that is fully based on the coordinate space HFB formalism, so that the theory can 
take into account coupling to the continuum configurations both in describing the pairing in the ground 
state and in description of collective excitations. A novel feature is that it includes, for the first time in 
the quasiparticle linear response formalism, the configurations where excited two quasiparticles are both 
occupying the continuum orbits above the threshold. We pay a special attention to consistency between 
the treatment of the pairing correlation in the static HFB and that in the linear response equation for 
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the collective excitations. This selfconsistency is quite important as we shall demonstrate in the following. 
We also discuss some basic features of pairing effects on collective excitations, taking as an example the 
quadrupole response in neutron rich oxygen isotopes including drip-line nucleus 24 O. 

The sections §§2 and 3 are devoted to derivation of the basic equations of the linear response equations. 
Numerical results for oxygen isotopes are discussed in §4, and conclusions are drawn in §5. 



II. TIME-DEPENDENT HARTREE-FOCK-BOGOLIUBOV EQUATION IN COORDINATE 

REPRESENTATION AND LINEARIZATION 



A. TDHFB in coordinate space 



The linear response theory or the random phase approximation can be formulated generally as a small 
amplitude limit of the time-dependent Hartree-Fock (TDHF) equation. Now in order to include the pairing 
correlations, we shall base our formulation on the time-dependent Hartree-Fock-Bogoliubov (TDHFB) theory. 
The TDHFB equations is often formulated with use of the discrete shell model single-particle basis and the 
Thouless representation |l||3^|, which however are not convenient to treat the continuum states. Thus we shall 
start our formulation with writing the basic equations of TDHFB in the coordinate space representation. It is 
an extension of the coordinate space HFB of Ref. to time-dependent problems. We follow the notation of 
Ref. in many aspects, while there are some differences since we do not impose the time-reversal symmetry 
assumed in Ref. [Q. 

The ground state of the system |$o) in the HFB formalism is a vacuum of the Bogoliubov quasiparticles. 
We denote the annihilation and creation operators of the quasiparticles by {f3i, /3j}(i — 1, • • •), which satisfy 
the vacuum condition j3i |<&o) — 0. The nucleon operators ip(r<j) and ip^{ra), where a — ±1/2 denotes the 
spinor component, can be expanded by the quasiparticle operators as 

^(r<7) = ^^(ra)A + ^M/3/ = ^Vi J i(r ff )A--^i(r5)4, (la) 

i i 

^{ra) = U*(ra)0l + V % (ra)^ = ^^(ra)/?? - (lb) 

i i 

Here ipi^(ra) and tp2,i(j&) are the single-quasiparticle wave functions that satisfy the coordinate space HFB 
equation |7j (see also the later description). We do not write explicitly the isospin degrees of freedom for 
simplicity although in actual applications the pairing correlation is taken into account separately for neutrons 
and protons. In the present paper we use functions (pi t i(ra) = Ui(r<r) and ip2,i(ra) = Vi(ra) instead of 
Ui(ro) and Vi(ra) for the convenience of notation whereas Ui(ra) and VAra) correspond more directly to 
the U, V matrices in the HFB formalism in the discrete shell model basis ||]. The symbol ip(ra) for a spinor 
function Lp{ra) denotes tp(ra) = (— 2a)ip(r — a) = (—i<j y ip)(r<7). For the time-reversal convention, we employ 
T(p(ra) = ip*{ra) = {—ia y Lp*){ra). These conventions are the same as in Ref. 0. We also use a notation 
ipi(ra) = ip*(ra) for the time-reversed function. 

In a time dependent problem, the system is described by a TDHFB state vector \&(t)}. We assume that at 
the initial time to the nucleus is in the ground state, i.e. |$(t = to)) = |$o)- Time evolution of the TDHFB 
state vector |$(t)) is governed by the time-dependent Schrodinger equation 

mim)) = kt)m)) (2) 

with the time-dependent mean-field Hamiltonian h(t), which can be expressed generally as 
h(t)= J J drdr'^2h{r(j,r'a' ,t)^{ra)ip{r'a') 

aa' 

+ i J J drdr' ^ \h(ra, r'a', t)^ (ra)^ (r'a') + h* (ra, r'a', t)ip(r' 'a')ip{ra) j 

-A J drJ2^(r°)ip{r<j), (3) 
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where the second and the third terms represent the pair potential, and A is the chemical potential or the 
Fermi energy. There are symmetry properties h*(ra, r'a', t) = h(r'cr' , rer, t) and h(r<r, r'cr', t) = h(r'&' , ra, t), 
which stem from the hermiticity of h(t) and the anti-commutation relation of the nucleon operators. 

Using the unitary operator U(t) that describes the time evolution of by |$(i)) = U(t) |$o), the nu- 

cleon operators in the Heisenberg representation ip^rat) = W (t)ip^ (ra)U(t) and ip(rat) = W \t)i/j(ra)U (t) 
are introduced. They also have expansion in terms of the complete single-quasiparticle basis, given by 

^{rat) = <Pi,i(r<rt)0i ~ </>2,i(™*)/??, ( 4a ) 

i 

divert) = <i( r ^l ~ V2,<(r*t)A. (4b) 

i 

The single-quasiparticle wave functions Lpi.i(rat) and ^^{rot) are now time-dependent, and they are at the 
initial time t = to set to the ground state quasiparticle functions ipi y i(ra) and ip2,i(r<j). The field equation of 

motion ih-^.ip(rat) = ip(rat),hH_(t) with huit) = (t)h(t)U (t) leads to the time-dependent Schrodinger 

equation for the single-quasiparticle wave functions tpi^(rat) and ip2,i(rat), that is written as 



d r 
ih—<t>i(rot) = / dr'22H(ra,r'a',t)^{r'a't), 



(5) 



where we have used a 2 x 2 matrix representation defined by 



- Hi , , , , h{ro,r'a',t)- \5{r ~r')5 aa , h(ra,r'a',t) . 

1 h*(ra,r'a',t) -h*(ra, r'a', t) + XS(r - r')6 ac ' ' 



and 

Note also that the same time-dependent single-quasiparticle equation (||) holds for </>j(rcrf) defined by 



The bar implies the operation of [ ^ Q )' With these notations, Eq.([4|) is written in a compact form 



\ ^(rdt) 



= Y,Hrvt)Pi + <k(rcrt)(3\. (9) 



A conjugate pair is formed by <pi(rat) and <pi(rat). The orthonormality condition and the completeness for 
the single-quasiparticle states are expressed as 

JdrY, 41™^ (rat) = J dr ^^(r<rt)^(rtrt) - S ij: (10a) 

a a 

dr^ 4 (rat)fy (rat) = 0, (10b) 
YUrvWHr'a't) +^(r<rt)^(rVt) = S(r - r')<W (\\\ (10c) 

i ^ ' 

It is convenient to use the time-dependent density matrices in describing the evolution of the TDHFB state 
vector. We define the normal and the abnormal density (pair density) matrices by 

p(ra,r'a',t) = ^(r'a')ip(ra) = ($ | ifi(r'a't)i)(rat) |$ ) , (Ha) 

p(rcr, r'cr', t) = (<Z>(t)\ il){r'a')ij){ra) |$(t)) = ($ | ^(r'a't)tp(rat) |$ ) , (lib) 
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respectively. These definitions agree with those in Ref. if the state vector time-independent = |$o) 

and the HFB ground state is time-even T |$o) = l^o), although in this paper such assumption is not made. 
Also our definition of h(ra,r'a' ,t) and h(ra,r'a' ,t) reduces to those in Ref. |?J with the same conditions. 
The two density matrices are combined in a generalized density matrix 1Z as 



(^(t)\^(r'a')ip(ra) |$(t)) ($(t)| ip(r'5 J )i/>(ra) |$(t)) 
<$(t)| ^(r'a')^(ra) |$(t)) ($(t)| ^(r , or / )^t( r 5-) |$(t)) 

p(ra,r'a',t) p(ra,r'a',t) 
p*(ra,r'a',t) 5(r - r')<W - p*(ra,r'a',t) 



(12) 



The generalized density matrix can also be expressed in terms of the time-dependent single-quasiparticle 
wave functions as 

Tl(ra,r'a',t) = ^2^(rat)^(r'a't). (13) 



B. Linearization 



Let us consider a TDHFB state vector |$(t)) which fluctuates around the HFB ground state |$o) under 
perturbation of an external field. The external perturbation induces also fluctuation in the selfconsistent 
mean-field. We now write the time-dependent mean-field Hamiltonian as h(t) — ho + V(t) where ho denotes 
the static selfconsistent HFB mean-field Hamiltonian for the ground state while the fluctuating field V(t) 
contains both the external field V ext (t) and the induced field V lnd (t). With presence of the pair correlation, 
the fluctuating field V(t) is a generalized one-body operator including pair creation and annihilation, which 
can be expressed as 



V(i) = J J drdr'^2{v(ra,r'cr',t)^(ra)ilj(r'a') 



+^v(ra,r'a',t)ip^(ra)^(r'a') + ^v*(ra, rV, t)ij)(r' a')ip(r' a')\ , (14) 

where we assume V(t) is hermite. To describe response of the time-dependent single-quasiparticle wave 
function <f>i(rat) to the perturbation, we expand it up to the linear order as 

cf>i(rat) = e -^(*-*o)/ft ( 0j ( r(T ) + Sfairat)) . (15) 

Here <pi(ra) and Ei denote the quasiparticle functions and the quasiparticle excitation energy defined as a 
solution of the static HFB equation for the ground state, 

H (j>i{ra) = Eifaira), (16) 

where Ho is the 2x2 matrix representation of the static HFB mean- field Hamiltonian Hq. 

An useful tool in describing the linear response is the single-particle Green function, which in our case is 
the single-quasiparticle Green function defined for the static HFB mean-field Hamiltonian ha. Because of the 
pairing correlation, we need both the normal and abnormal Green functions 

G (rat, r'a't') = -i9(t - t') <$ | {^(rat), ^(r'a't')} |$ ) , (17a) 
F (rat, r'a't') = -iO{t - t') ($„| {^(rat), ^(r'a't 1 )} |$„) . (17b) 

Here the nucleon operators ^(rat) = U^t)^ (ra)Uo(t) and tp(rat) — Uo(t)ip(ra)Uo(t) with Uo(t) = 
e -i{t-t )ho/n evo i ve m time under the static HFB mean-field Hamiltonian ho. More useful is the single- 
quasiparticle Green function in the 2x2 matrix representation |35|,[36| defined by 

G (rat r'a't'\--if)(t-t') ( (®o\ H(rat),^ (r'a't')} \$ ) (<S> \{i,(rat),iP(r'a't')}\$ ) \ (m) 
y o[ rat,rat)- iU[t t ) I ^ {^i ( r at), ip^ (r'a't')} |$ ) ($o| {ip f (rat),ip(r'a't')} |$ ) / ' [ ' 
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The normal and abnormal Green functions Go and Fo are contained in the 11 and 21 components, respectively. 
In the present paper we use the retarded functions in stead of the causal or Feynman functions. We call Qq 
the HFB Green function. The Fourier transform of the HFB Green function Qo(rat,r'a't') is given by 
Qo(ra, r'a', E + ie) = (E + ie — Ho) where the infinitesimal constant e — > +0 exhibits the causality. In 
actual numerical calculations, e is fixed to a small but finite number. Note that the HFB Green function can 
be constructed in different ways. If we adopt the spectral representation, it is expressed as 

5o( „, rV , E ) = S M|)^ + «!|)3|V). (19) 

i 

Another form which is useful to treat the continuum states is given in Ref. which we utilize in the 
following. 

The Fourier transform of the linear response 5<pi of the single-quasiparticle wave function is now expressed 

as 

5(f> i (ra,uj) = J J dr'dr" Go(ra,r'a' Juj + ie + E i )V(r'a',r"a",u})(f> i (r"a") (20) 

a' a" 

in terms of the HFB Green function. Here V is the matrix representation of the fluctuating field V(t), which 
is defined in the frequency domain by 

V(ra, r'a', u ) = ( <J°2 "\ ^'/^ s ) ■ (21) 
v ' ' ' \v*(ra,r'a',u>) -v (ra, r'a , w) J v ; 

For the conjugate wave function (fa holds the same equation except that the energy argument of the HFB 
Green function is replaced by Tilj + ie — Ef. 

Using the above results, the linear response in the density matrix lZ(u>) = + 51Z(u>) is obtained as 



6K(ra,r'a',uj) =11 dr 1 dr 2 ^ ^Qo(ra,riai, -E{ + Tiuj + ie)V(riai,r 2 a 2 ,w)<j)i(r 2 a 2 )(j)\{r' a') 

+4>i{ra)4>\{r 1 a 1 )V(r 1 ai,r 2 a 2 , uj)Q (r 2 a 2 ,r' a' , ~E { -huj- ie) j . (22) 

Expectation value A(t) = (&(t)\ A \<&(t)) of a one-body operator 

A = i /" J drdr'^2 {An(ra,r'a')^(ra)i/j(r'a') + A 22 (ra,r'a')i/j(ra)^(r'a') 

+A 12 (ra,r'a')^ (ra)^ (r'a') + A 21 (ra,r'a')^(r'a')i(;{r'a')} (23) 
is expressed in terms of the generalized density matrix as 

A(t) = J f drdr'^2TrA(ra,r'a')n(r'a',ra,t) (24) 



with 



A(ra,r'a>) = ( A /f a ' r '/,\ tf"' ) • ( 25 ) 

v ; \A 2 i(ra,r'a') A 22 (ra,r'a') J' y ' 

where Tr denotes the trace with respect to the 2x2 matrix. Linear response of A(t) in the Fourier represen- 
tation is given by 

5A(ui)= J " J dridr 2 dr 3 dr^ ^ ^ 

|Tr^l(riCTi, r 2 a 2 )g Q (r 2 a 2 , r 3 a 3 , -Ei + Tiuj + ie)V(r 3 a 3 , r 4 cr 4 , Ci/)<^(r 4 o- 4 )4 (rid) 

+TrA(riai,r 2 a 2 )(/)j(r 2 a 2 )(l)l(r 3 a 3 )V(r 3 a 3 ,r4,a4 : ,LL!)go(r4 : a4 1 riai,-E l - %u> - ie)} . (26) 
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C. Unperturbed response functions for local densities 

In the following we assume that the fluctuating field V(t) including both the external and the induced 
fields is expressed in terms of the spin-independent local density 

p(r) = ^ (ra)ifj{ra) (27) 

and the local pair densities 

PHr)=±J2^( riT )^( r °)> ( 28a ) 

P(r)=^J2i)(ra)iP(ra). (28b) 

a 

The above assumption is valid when the residual nuclear interaction is a zero-range force. For the convenience 
of notation, we use a symmetrized form for the pair densities, defined by 

p + (r)= P f (r) + P(r), (29a) 

p_(r)=pi(r)-P(r). (29b) 

The fluctuating field is expressed as 

V(t) = Jdr (wo(r,t)p(r) + w + (r,t)p + (r) + «_(r,t)p_(r)) . (30) 
We consider now linear response in expectation values of the local density operators 

p(r,t) = ($(t)|#r)|$(t)> = X>(r<r,rcr,t), ( 31a ) 

(7 

p + (r, t)= ($(t) | p + (r) |*(t)) = \ (p* (ra, ra, t) + p(ra, ra, f)) , (31b) 

(X 

P-(r,t)= !*(*)) = \Y,{f>*{™,ra,t) - p{ra,ra,t)) , (31c) 



under the influence of the fluctuating field V(t) and the static HFB mean-field Hamiltonian ho. The relation 
between V{t) and the density linear responses is expressed as 



Sp a (r,t) = f dt> [ dr'Y,^{rtyt')v {r',t') 



(32) 



by means of the unperturbed response function Rq (rt,r't'), whose formal expression is p(| 

Rf(rt,r't') = -i9(t-t') (*„| [u^t) Pa (r)U (t),U^t')pp(r')U (t')] |$ ) . (33) 

Here p a (r,t) represents the three kinds of density p(r,t), p+(r,t), and p_(r,t). We also denote p a (r) = 
p(r), p + (r), p_(r), and v a (r,t) = v (r,t),v + (r,t),v-(r,t), respectively. With use of Eq.(p6|), the Fourier 
transform of the unperturbed response function is obtained as 

Rf(r,r',u) = ^YIY, {^H^ofw, r'a', -E, +hcj + ie)B&(r'a') 

i aa' 

+^l(r'a')BG (r'a / ,ra, -E, - huj - ie)A^(ra)\ (34) 



G 



in terms of the HFB Green function Qq. Here A and B are 2x2 matrices which correspond to the normal and 

o o ) ' ( o i ) ' and ( 1 



abnormal densities; A = ( 2 n ° n ) , ( I ? ), and ( ? i ] for Pa (r) = p(r),p + (r), and p_(r), respectively. 



B = ( q ^ J for pp(r) = /5(r), and S is the same as ^4 for /5^(r) = p + (r), p_(r). 

If we insert Eq.()l9|) to the above equation, the unperturbed response function reduces to a spectral repre- 
sentation in a standard form 



Rf(r,r',u;) = ~ £ |(0| p a (r) \ij) (ij\p (r') |0) 



-<0|^(r')|ij) (ii|Pa(r)|0) 



where 



(ii|p a (r)|0)=^4(ra)A^(ra), 

(T 

(0|p a (r) \ij) =^^H^(m). 



1 



+ ie — Ei — Ej 
1 

fiu + ie + E, + E, 



(35) 

(36a) 
(36b) 



This can also be derived directly from Eg. (|33|) . This expression, however, is not convenient to treat the 
continuum states. 



III. CORRELATED LINEAR RESPONSE WITH CONTINUUM STATES 

A. Integral representation 

If the nucleus is put in the space of infinite volume, the spectrum of HFB single-quasiparticle states become 
continuum for the quasiparticle excitation energy E which exceeds the one-particle separation energy Si = |A| 
0,^]. Accordingly, the single-quasiparticle HFB Green function Qq also exhibits the continuum spectrum. 
Furthermore, the HFB Green function at the continuum energy is required to satisfy the boundary condition 
of outgoing wave in the exterior region of the nucleus |3S| ], which is a suitable boundary condition for the 
continuum states. Let's us implement these features in the response functions. 



(a) 
ie' 
-i e' 

(b) 



C 



-x — x 



ImE 



ReE 



-x — x- 



-Ecut 
C 







-X 



X 



-x — x- 



-x — x- 



FIG. 1. Contours C and C in the integral representation of the response function. The crosses represent the poles 
at E — ±Ei corresponding to the bound quasiparticle states. The thick lines are the branch cuts. 

This task is accomplished in two steps. Firstly, we handle the continuity of the HFB quasiparticle spectrum. 
For this purpose, we rewrite the summation over the HFB quasiparticles in Eq.(|34|) with use of an integral 

representation. Namely, a summation /(— £"i)^(r)^ (tct) is replaced by a contour integral in the complex 
E plane as 



** c 



(37) 
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in terms of the HFB Green function Go, Eq.([l9|), which has poles at E — — Ei on the negative energy axis. 
The contour has to be chosen so that it encloses these poles, and at the same time, it must avoid the poles 
associated with Go(E + Hlu + ie) and Go(E — %u> — ie) that account for Go(—Ei + hu + ie) and Go(~Ei — %uj — ie) 
in Eq. (|34|) . A contour C satisfying this requirement is shown in Fig. 1(a), where the constant e' must satisfy 
e > e' > 0. Consequently, the unperturbed response function is expressed as 

Rf{r,r',u>) = ^ ( dEy"{TrAg (ra,r'a',E + hLJ + ie)BGo(r'a',ra,E) 
47ri Jc , 

(J (J 

+TrAGo{ra, r'a', E)BGo(r'a', ra, E-hu- ie)} . (38) 

The contour is the same as those adopted in Refs. ]35|,[l7]] except the condition for e' . 
For spherically symmetric systems, the partial wave expansion can be applied as 

0(nr) = y yro (m)^ y (r), 0y(r) = ( ) , (39) 

G (ra,r'a',E) = ]T Y ljm (ra)^Go,lj (r, r', E)l£ m (rV), (40) 

Jif (r,r» =E F ^^)^72^o,i(nr',o;)y i V(P), (41) 

with Yy m (fcr) being the spin spherical harmonics. The unperturbed response function with the multipole L 
is then given by 



i$£(ry,«) = -L / ^^o FT {Tr^go,^-- (r, /, E + ^ + ie)BGo, lj {r\ r, E) 

+TtJGo,ij(r, r', E)BGo,v?{r' >, E-hw- ie)} . (42) 



If we treat separately the discrete and the continuum parts of the HFB spectrum, the following equivalent 
expression is obtained 

iZg (r, r', W ) = ~ {V3 'lI L +i )2 ffii M- 4 ^' r '' + ^ + i£ )^«« ^ 

+4>lij{ r ')BGo,i>j>{r' ,r,-E nlj - fiu; - ie)A~$ nlj (r)| 

+TtAGo,i j (r, r', E)BGo, Vj > (r', r,E-hu-ie)}, (43) 



where the summation ~Y^' n runs only over the bound single-quasiparticle states with discrete spectrum E n ij < 
|A|, while the contour C" in the complex plane encloses only the continuum part, as shown in Fig. 1(b). 

Having the integral representation of the response functions, we impose on the single-quasiparticle Green 
function the boundary condition of out-going wave in the exterior region. The HFB Green function satisfying 
this requirement is given in Ref. |35|Jl^ ], which can be applied to spherically symmetric systems with local 
potential. Namely, the HFB Green function is constructed |35|Jl7|] as 

Gom (r, r', £0 = E <tf W " ^ ^If ° V, E) + 6(r> - r)4f\r, E)^ s)J \r> , E)) . (44) 

s,s'=l,2 

Here 0y (r, E)(s = 1,2) are two independent solutions, regular at the origin r — 0, of the radial HFB 
equation 
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while (f>ij{r, E)(s = 1, 2) are two independent solutions that satisfy the boundary condition at r — > oo 

(r, E) -> ( ^ ) , (r, £0 - ( eife _% )r ) • (46) 

Here k±(E) = yj2m(\ ± E)/T% and the branch cuts are chosen so that Imfc± > is satisfied. Thus the 
nucleons at r — > oo is out-going for E + ie at the continuum energy E > |A| [ ^5|Jl7| . The coefficients of? (E) 
are expressed in terms of the Wronskians as 

with 

(48) 



B. RPA response functions 

Let us now describe linear response of the system by including correlation effects caused by the residual 
interactions among nucleons. The correlation brings about the induced field, which can be described as a 
fluctuating part of the selfconsistent mean-field associated with the TDHFB state vector |$(t)). Generally, 
the expectation value of the total energy for |$(t)) is a functional E[1Z] = E[p, p, p*] of the generalized density 
matrix 1Z, or p, p, and p* . The selfconsistent mean-fields h and h in Eq. (^|) are then defined by a derivative 
of the functional with respect to the densities, 

OE 

h(ra, rV, t) = —— (49a) 

opyr'a' , ra, t) 

r)E 

h(ra, rV, t) = 2— — -. (49b) 

op*[r'a, ra, t) 

If we assume zero-range effective interactions, the energy functional is expressed in terms of the local spin- 
independent densities p a (r,t) — p(r, t),p + (r, t) and p-(r,t). In this case, the selfconsistent mean-field is 



expressed in the same form as Eq.(30) with the field functions given by 



, \ 9E 

Mr,t) = 50 
op a (r,t) 



Accordingly, the induced fields are 



Inserting v a = v" ld + v e * 1 to Eq.(|32]), we obtain the equation for the density linear response 5p a (r,t) = 
6p(r,t),6p+(r,t), and 5p-(r,t), with the correlation effect taken into account. That is, 



Sp a (r,u>) 



J dr£<V,r» |^V)^ 7 (r>) + ] ■ (52) 
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This linear response equation is similar to the one in the continuum linear response theory for unpaired 
systems [^J, but we here take into account the fluctuations in the pair densities. It is possible to write the 
above equation as 

3u 

8 P {uj) = R (lu)—5 P (lu) + R Q (uj)v ext , (53) 

where the three kinds of density fluctuations are represented as a single extended vector Sp(uj) — 
(5p(r,uj),Sp+(r,uj) 1 5p-(r,uj)) T with three components, and the same representation applies to the external 
fields v ext . We consider both neutrons and protons in the actual application. Then the density response dp 
have six components (three for each isospin). The solution of the linear response equation is expressed as 

Sp{u) = (i-Ro(lj)^] R a (uj)v ext = R(cu)v ext (54) 



dp 



or 



Sp a (r,uj)= / dr'Y l K xf> {rJ>u)vT{r% (55) 

J 

where R(uj) = [R a ^(r, r', to)) is the correlated response function. This is nothing but the response function 
in the random phase approximation (RPA). 

The density response of the spherical system to an external field with the multipolarity L is given by 

v e *\r)=Y LM {r)v aL {r), (56) 

Sp a (r,u) = Y LM (r)-^5p a L(r,uj), (57) 

and the linear response equation 

Sp aL (r,uj) = J^dr'R%f L (r,r',w) (j^(r')-^ 5p aL (r', u) + v aL (r' , uj)^ . (58) 

This equation can be solved numerically with use of the mesh representation of the radial coordinate in a 
way similar to Ref. || . The strength function for the response to the external field is expressed by means of 
the RPA response function as 

k 

= _Il m f f drdr'J2R%(r,r\u>)v* aL (r)vp L (r') = --Im f dr V < L (r)<5p ai (r, w ). (59) 

Let us discuss characteristic features of our linear response theory by comparing with other approaches. 
Firstly, the present formalism includes the continuum quasiparticle states wherever they appear. Previous 
continuum QRPA's (or the continuum quasiparticle linear response theories) p6|-|29]] adopt approximations 
that take into account the continuum boundary condition only for the normal Green function Go (E) , but not 
for the abnormal Green function Fq(E), while we have treated the single-quasiparticle HFB Green function 
exactly by using the construction by Belyaev et al. J3f| . We also emphasize that the present theory includes 
the configurations where two nucleons are both in the continuum states. Since the HFB single-quasiparticle 
spectrum is made of the discrete bound orbits and the continuum unbound states, two-quasiparticle states 
are classified in three groups. The first is the one where both quasiparticles occupy two discrete bound 
states. The second is where one quasiparticle is in a bound state with discrete excitation energy Ei while 
the second quasiparticle occupy the continuum unbound states whose excitation energy E exceeds |A|. These 
one-particle continuum states emerge above a threshold excitation energy E t h.i = minEi + |A|, where minEi 
is the minimum quasiparticle excitation energy. In the unpaired limit, this threshold energy becomes identical 
to the excitation energy —eh,iast of nucleon from the last occupied hole orbit to the zero-energy threshold. The 
last group of two-quasiparticle states is the one where both two quasiparticles occupy continuum unbound 
states. The threshold excitation energy for the two-particle continuum states is E t h,i = 2 j A j , i.e. it is 
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twice the one-particle separation energy Si- The present theory includes all the three groups as seen in 
Eq. (|43|), whereas the previous approaches include only the configurations where only one nucleon is in the 
continuum states. Escaping processes of one particle to the external region takes part in at the excitation 
energy Tim > E t h,i above the one-particle threshold, and consequently the strength functions show continuum 
spectra. At excitation energy hu> > E t h,2 above the two-particle threshold, processes of two-particle escaping 
contribute. A related point is that the previous continuum QRPA's adopt the BCS approximation that has 
difficulty to describe the pairing in the continuum, whereas we have exploited the coordinate space HFB 
formalism to remove the shortcoming. 

Secondly, the present theory takes into account the particle-particle correlations associated with fluctua- 
tions in the pair densities <5p_|_, and <5p_, which always emerge in paired systems. In other words, the pairing 
residual interaction contributes to the linear response equation (p8|). This particle-particle correlation has 
been neglected in many continuum quasiparticle linear response theories or QRPA's. Even when it is in- 
cluded |2f|^7||, its effect is not clarified. We emphasize that this contribution should be included to keep a 
selfconsistency of the pairing correlation, and it indeed has important consequences. This point is discussed 
in detail in the next section. 

Thirdly, the present theory reduces to the continuum linear response theory for unpaired systems || if 
the HFB ground state has the zero pairing potential. This is easily seen in Eq. ([f3|) , which, in the zero 
pairing limit, has only the discrete part consisting of the occupied hole orbits 4» n ij- Note also that only the 
11 component (the normal Green function Go) of Go contributes to the response function for the density 
operator p(r). In this sense, the present theory is an extension of Ref. || to paired systems described by the 
coordinate space HFB. 

IV. NUMERICAL ANALYSIS 

We apply the above formalism to spin independent quadrupole excitations in open-shell oxygen isotopes 
near the neutron drip-line, where both the pairing correlation and the continuum effects may play important 
roles. In the following, we mainly discuss results for the drip-line nucleus 24 0. Since our purpose in this 
paper is to clarify basic characteristics of the theory, rather than to make a precise prediction, we adopt a 
simple model Hamiltonian which consists of a spherical Woods-Saxon potential and the residual two-body 
interactions. 

As a pairing force we assume the density-dependent delta force [pl,j37| 

v pair (r,r') = lv (l-P a )(l-p(r)/p )S(r-r'), (60) 

where P a is the spin exchange operator. The interaction strength depends on the position through the total 
nucleon density p(r) = p n (r) + p p (r). With po = 0.16fm~ 3 (the saturation density), the pairing force is more 
effective at low-density nuclear surface than in the interior region. The density-dependent pairing force is 
widely adopted in the HFB calculations not only for unstable isotopes P, p^ , p^| , ^9| , p2p^ ] but also for rapidly 
rotating deformed nuclei (|^^] . 

To obtain the ground state, we solve the radial HFB equation with the box boundary condition cf>i(r = 
Rmax) — according to the procedure of Ref. Q], and determine selfconsistently the nucleon density p q (r) (q = 
n,p), the pair density p q (r) and the pairing potential A g (r) = Vb(l — p(r)/ po)p q (r). For the Woods-Saxon 
potential we adopt the standard parameter set ^9|j40t| . We have used the radial mesh size Ar = 0.2 fm, 
and the cut-off at R max — 20 fm. To evaluate the densities and the pairing potential, we include all 
the quasiparticle states whose excitation energy is below E max — 50 MeV and with the orbital angular 
momentum l max < 7h. The pairing force strength Vq is chosen Vo.„ = Vo tP = 520fm~ 3 MeV so that the average 
neutron pairing gap (A n ) in 18 ^ 20 O approximately agrees with the overall systematics A w 12/\/A ~ 2.8 
MeV Here we define the average gap by (A„) = J drp(r)A n (r)/ J drp(r), which corresponds to 

(A n ) uv — ^2 u^v^A^/ J2fj, u n v n (a 1 denoting the canonical basis) adopted recently in the literature pl| , |38| , p2| . 
The calculated average neutron gap is (A„) = 2.74,3.13,3.30, and 3.39 MeV for i8,20,22,24Q ; respec ti V ely. 
Note that the pairing gap increases as approaching the neutron drip-line line. For comparison, we use 
also the density- independent pairing interaction (the volume- type pairing), which is given by dropping the 
density-dependent term in Eq.(|60|). In this case, we use V = 240fnV 3 MeV, po = oo so that the calculated 
average neutron gap (A n ) = 2.70,2.69,2.26, 1.64 MeV (for 18 > 20 > 22 > 24 O) gives approximately the same value 
in 18,20 O. The neutron pairing gap for the volume-type pairing decreases as increasing the neutron number, 
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especially at the drip-line nucleus O. This trend is opposite to that for the density-dependent pairing, 
and consequently the calculated pairing correlation differs significantly for 22 > 24 located near drip-lines even 
though the force parameter is adjusted in more stable 18,20 O. Apparently, the density-dependent pairing force 
favors the pairing correlation in the low-density surface region which develops in nuclei near the drip-line. 
This is a characteristic feature of the density-dependent pairing |?],|],[l4j ■ 

The quadrupole response of the system is described by means of the linear response formalism presented 
in the previous section. As the residual interaction acting in the particle-particle channel, we use the same 
pairing interaction as that used for the static HFB calculation in order to keep the selfconsistency, which 
turns out very important in the following. For the particle-hole channel, we assume the Skyrme-type delta 
force 

v ph (r,r') = {t {l + x P a ) + t 3 (l + x 3 P a )p{r)) 5{r ~ r'), (61) 

where we adopt the same parameter as Ref. @; t = / x (-1100)fm 3 MeV, t 3 = f x 16000fm 6 MeV, x = 0.5 
and X3 = 1. Here the renormalization factor / is adjusted by imposing a selfconsistency Q that the calculated 
static polarizability for the isoscalar dipole field D[fj = J2 n r ^iM + J2 P r ^iM becomes zero for each nucleus. 
The linear response equations are solved in the following way. The unperturbed response function is obtained 
by evaluating Eq.(^2|) where the contour integral is performed numerically with an energy step \AE\ = e/16. 
We use a small smoothing constant e = 0.2 McV, which corresponds to a smoothing with a Lorentzian 
function with FWHM of 0.4 MeV. As for the HFB Green function we adopt the exact expression (fl4"|), and 
evaluate the radial functions (also the response functions) by using the mesh Ar = 0.2fm and the cut-off 
radius R max — 20fm. In evaluating Eq. (ft2|) , we introduce also the angular momentum cut-off l ma x — 7h. 
These parameters are the same as in solving the static HFB calculation. The cut-off energy for the contour 
integral is E cut = 50 MeV, which is a natural choice since it is the same cut-off energy E max for the static 
HFB. Because the value of E cut is larger than the potential depth (measured from the Fermi energy), the 
contour integral can include all the hole orbits in the limiting case of the zero pairing potential. We thus 
guarantee that in the zero pairing limit the evaluated response function for the density operator becomes 
identical to the response function for unpaired systems [B. The RPA response functions R a ^(r,r' ,huj) is 
obtained by solving Eq.([58|), which is a 600 x 600 linear matrix equation in the radial mesh representation. 

Figure ^ shows the strength functions calculated for the quadrupole excitation in 24 0. Here we show the 
strength function S T (huj) = dB(Q T 2,0^ s — > k)/d{Tiu) = J^Mk (0gs| Q2M \k) 2 S(Tiu> — frujk) for the neutron 
and proton quadrupole moments Q% M — J2n r2 ^ 2M an< ^ Q2.Ni ~ J2 P r2 ^2M, the isoscalar and isovector 
quadrupole moments Q ! 2 m — Q2M + ano - Q2X4 = Q2M ~ Q2M as a function of the excitation energy Hiv. 

Plotted here is the strength functions summed over the magnetic quantum number M so that their energy 
integral are equivalent to B(Q T 2,0+ S — > 2 + ). A noticeable feature seen in Fig.|^ is presence of a sharp and 
intense peak around hui « 5.0 MeV. There also exist a resonance peak around Tiui « 17MeV with a width of 
~ 3MeV, and in addition a very broad distribution around huo sa 20 — 40 MeV. These two modes correspond 
to isoscalar(IS) and isovector(IV) giant quadrupole resonances (GQR) whereas they are not pure ISGQR nor 
IVGQR since a large admixture of isovector (isoscalar) strength is present in each resonance region of the 
isoscalar (isovector) modes. 
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FIG. 2. The strength functions for the quadrupole moments QJm ( r — IS,IV,n,p) in 24 with use of the 
density-dependent pairing force. The renormalization factor for the particle-hole residual interaction is / = 0.632. 
Solid and dotted curves represent those for Qifi, Q2M) and Q%m > Qsjmi respectively. The inset shows a magnified part 
for the low excitation energy hui < lOMeV. The threshold energy E t h,i = E t h,2 = 2|A| = 4.14 MeV is indicated by an 



The peak excitation energy of the low-lying mode is Tiuj — 5.0MeV. Since the threshold energies for one- 
neutron and two-neutron continuum states are, E t h,i — Eth,2 — 4.14MeV, the low-lying neutron quadruple 
mode is embedded in the neutron continuum states. The FWHM width of this mode w 45keV is almost the 
same as the smoothing width 2e = 40keV, indicating that the mode has very small escaping width in spite 
of the coupling to the neutron continuum states. 

The low- lying quadrupole mode at Tilo « 5.0MeV is characterized by very large neutron strength. The 
isoscalar strength is also large, but the proton strength is small. The isovector strength is sizable but smaller 
than the neutron and isoscalar strengths because the small proton contribution is in phase with the neutrons. 
The integrated neutron strength B(Q n 2) below Hlu = 6.0 MeV is 605fm 4 whereas B(E2)(= e 2 B(Q p 2)) is 
only 16.4e 2 fm 4 . Converting the strengths to the neutron proton ratio M n /M p for the quadrupole transition 
amplitudes, it is evaluated as M n /M p =6.1. The neutron character of this mode is apparent by comparing this 
value to a simple macroscopic estimate M n /M p = N/Z = 2. The low-lying mode exhausts approximately 
20%, 10% and 19%, of the sum rule value for the isoscalar, isovector, and neutron quadrupole transition 
strength, respectively (See Eq.(|62|) and Fig. [fa]). This ratio is much larger than a typical value for the 
low-lying isoscalar quadrupole mode (~ 10%) Q. On the contrary, the proton strength (i.e. the B(E2) 
strength) is small, exhausting only 2% of the energy weighted sum rule. The strength of the low-lying neutron 
mode is continuated to the neutron strength distribution in the interval (Hlu « 6 — 15 MeV) between the 
low-lying mode and the giant resonances, which corresponds to the so called threshold strength. However, 
the character of the low-lying neutron mode clearly differs from that of the threshold strength [g] since its 
large neutron strength of collective nature is generated by the residual interaction (see below). We note also 
that the HFB single-quasiparticle states have no discrete bound orbits, instead they are all continuum states 
or resonances. The low-lying neutron mode is a collective state made of continuum two-quasiparticle states. 

We have also calculated the quadruple strength functions for 18 > 2 °. 22 0. The peak energies of the low- 
lying mode is hui = 4.2,4.1,4.5,5.0 MeV for A = 18,20,22,24, respectively, which is 1 - 2 MeV higher 
than the experiments |32fl4|Q. The calculated neutron strength {B(Q n 2) = 152, 295, 426, 605fm 4 for 
A=18-24) increases sharply with increasing the neutron number. On the contrary, the calculated B(E2) = 
15, 18, 17, 17e 2 fm 4 (yl = 18 — 24) stays almost constant. The calculation underestimates the experimental 
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B(E2) value Uglgg by a factor of 3-1.5 for 18 ' 20 O, but with reasonable agreement in 22 0. The calcu- 
lated B{E2) values differs by a factor of 0.5-1.5 from the quasiparticle RPA calculation in Ref. |^l| , ^| . The 
adopted model Hamiltonian, especially the Woods-Saxon potential, should be improved to make more quan- 
titative comparison, e.g. by using more realistic Skyrme interaction both for the Hartree-Fock potential and 
for the particle-hole interactions. Leaving such improvements for a future investigation, we in the following 
focus on basic aspects of the present theory. 



300 



250 — 



-T50O 



> 

Si 200 



£i 150 
g 100 
50 








10 



20 



30 



40 



50 



fico [MeV] 



FIG. 3. The strength functions for the isoscalar and proton quadrupole moments Q21I1 and Q\ M in O with use of 



the density-dependent pairing force, calculated with the discretized continuum states (see text) . The renormalization 
factor for the particle-hole residual interaction is / = 0.702 
QIm> respectively. 



Solid and dotted curves represent those for Q2M an d 



The present theory takes into account the one- and two-particle continuum states without making any 
discretization method or bound state approximation. It is interesting in this respect to see difference between 
the results presented above and results that would have been obtained by using discretized continuum states. 
We show in FigJH the strength function obtained by using the discrete spectral representation Eq. instead 
of Eq.(^|). Here we solve discrete single-quasiparticle states with the box boundary condition (j>i{R m ax) = 0, 
and adopt all the quasiparticle states with the excitation energy Ei below the energy cut-off E max = 50 
MeV and the angular momentum cut-off l m ax = 7, including the discretized states in the continuum region 
Ei > |A|. The renormalization factor / is slightly changed in this approximate calculation in order to 
reproduce the zero energy dipole mode. With use of the discretized continuum approximation, the strength 
function consists only of discrete states (but smeared with the smoothing width) even in the continuum region 
above the threshold. This causes spurious fluctuation, as seen in Fig.0. It is noted on the other hand that a 
gross profile of the strength distribution is described fairly well by the discretizing approximation. It can be 
expected that in the limit of large box size R ma x — * 00, the discretizing approximation would give the same 
results with the one obtained with the continuum linear response theory, whereas such limiting is practically 
very difficult. Note also that the strengths associated with the low- lying neutron mode at hoj w 5 MeV are 
almost reproduced by the discretizing approximation. This is because the low-lying mode in the present 
model calculation is a very narrow resonance, for which a bound state approximation can be a reasonable 
approximation. 
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FIG. 4. Dependence of the quadrupole strength function in 24 on the pairing interaction. Plotted are the 
strength functions for the isoscalar quadrupole moment Qifi for the density-dependent pairing (dashed curve), for 
the volume- type pairing (solid), and for the case where the pairing is neglected (dotted). See text for details. The 
rcnormalization factor for the particle-hole residual interaction is / = 0.702. 

The pairing correlation play various important roles for the response. We show in Fig.^ the calculation 
without the pairing correlation (i.e. the pairing interaction is neglected both in the static HFB and the linear 
response equation). In this case, the quadrupole strength in the low-energy region becomes much smaller. 
The small peak at % m 3.0 MeV is not very collective, but it is basically non-collective neutron particle-hole 
transition 2s^ — > l<i| with slight enhancement due to correlation. It is clearly seen that the pairing increases 
the collectivity of the low-lying neutron mode. The pairing effect on the low-lying isoscalar quadrupole 
mode in stable nuclei is well known |j42| . Our result suggests a similar effect on the neutron mode, which is 
embedded in the continuum states in the case of drip-line nuclei. 

We found further that the low-lying neutron mode is quite sensitive to the density-dependence of the 
pairing. Figure ^ shows also the result that is obtained with the density-independent volume-type pairing 
force in place of the density-dependent one Eq.(|60|). The neutron (and isoscalar) strength in the low- lying 
neutron mode obtained in this calculation is significantly small, i.e. it is about 50% of that with the density- 
dependent pairing. Using this sensitivity, one may be able to probe the density-dependence of the pairing. 
The difference between the density-dependent and the volume pairing is most significant in 24 0, but not very 
large in the other oxygen isotopes, which is in accordance with the behavior of the average neutron gap (A„) 
discussed above. Fig. ^ shows that the pairing correlation influences also the strength distribution in the 
giant resonance region. It has an effect to reduce the strength in the high lying modes. An apparent reason is 
that the low-lying neutron mode collects nearly 20% of the total energy weighted isoscalar and neutron sums 
in the case of the surface pairing, and removes corresponding amount of strength from the giant resonance 
region because of the sum rule. The peak position and width (distribution profile) of GR's are also slightly 
affected by the pairing. 
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FIG. 5. The same as FigM but the dynamical pairing correlation effect in the linear response equation is neglected. 
Namely the density-dependent pairing is taken into account only in the static HFB mean-field. 

There are two kinds of mechanisms in the pairing effects on the collective excitations. The first is a 
static effect. Since the static pair potential A(r) in the HFB mean-field influences strongly the quasiparticle 
excitations, it hence affects the collective excitations. The second is dynamical. It is a correlation effect 
that is caused by the residual pairing interaction entering in the linear response equation (p8[). This arises 
because the collective excitations in the paired system induce not only the fluctuation in the density bp 
but also the fluctuations in the pair densities <5/5+ and 5p~, for which correlation is brought by the residual 
pairing interaction (see, Eq.([58|)). To quantify this dynamical pairing correlation effect, we show in FigJ^ a 
calculation in which the pairing interaction is neglected in the linear response equation (|58|), but included 
in the static HFB mean-field. Comparison with Fig.^ shows that the dynamical pairing correlation has a 
sizable effect to lower the excitation energy of the low-lying neutron mode by about lMeV. (One may also 
note that the width of the low lying mode increases slightly in this calculation. ) Performing a calculation 
where the particlc-holc interaction is neglected, wc find that the low-lying neutron mode is produced also by 
the pairing interaction alone, although the particle-hole residual interaction together is essential to make the 
collectivity large. If we used the the conventional monopole pairing force (the seniority force) instead of the 
zero-range pairing force, the dynamical pairing correlation effect would be missed. The dynamical effect is 
rather related to effects of the quadrupole pairing that has been discussed in connection with the low-lying 
quadrupole mode in stable nuclei. 
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FIG. 6. (a)Energy weighted sum J " hiv' S ISJV ' n ' p (hu}')hdio' of the quadrupole strength functions in 24 0, shown 
in Fig.^. Here the density-dependent pairing is taken into account selfconsistently both in the static HFB mean-field 
and in the linear response equation. The dashed horizontal lines show the sum rule values for isoscalar (the same for 
isovector) transitions, and those for the neutron and the proton transition strengths, (b) The same as (a), but the 
dynamical pairing effect is neglected (cf. Fig.^). Namely the pairing is taken into account only in the static HFB 
mean-field. 



One of the most important characteristics of the present theory is that the energy weighted sum rule 

f °° hwS^i^hdLJ = g (N <r 2 ) n + Z (r>) p ) , (62a) 
/" huS n -P(huj)hduj = {r 2 ) n , ^Z{r\, (62b) 

is satisfied very accurately. The sum rule should hold since the adopted residual interactions (both the 
particle-particle pairing force v pa i r and the particle- hole interaction v p h) keep the Galilei invariance ^2j, i.e., 
they commute with the quadrupole operator t 2 Y2m(j : )- This is demonstrated in Fig. ^(a) that shows the 
calculated energy weighted sum / tku'S T (hLo')hduj' (r = IS, IV, n,p) as a function of the energy boundary 
hu. We stress here that the selfconsistent treatment of the pairing both in the static HFB and in the 
dynamical linear response is crucial to satisfy the sum rule. This is shown in Fig. |(b). This plots the energy 
weighted sum obtained in a truncated calculation in which the residual pairing interaction is neglected in the 
linear response equations (cf. the corresponding strength function is shown in FigJE]). In this calculation, the 



sum rule Eq.(62) is strongly violated by about 50% for neutrons. The violation is significant especially in 
the high lying region (E^Z 10 MeV), as the strength in this region (shown in Fig.||) apparently overestimates 
that of the full calculation (Fig.^J) . The violation originates from the inconsistent treatment of the pairing 
interaction in the calculation of Figs.[| and ^(b), where the pairing correlation is included only in the HFB 
static mean-field, but neglected in the linear response equation. One can understand this by noting that the 
static pair potential A(r) violates the Galilei invariance, i.e., it does not commute with the density operator. 
By taking into account selfconsistently the pairing interaction in the linear response equation, the Galilei 
invariance of the original pairing interaction is recovered, and the sum-rule is satisfied. This corresponds 
to a previous result p7| showing that the sum rule violated by the monopole pairing potential can be 
remedied by including the selfconsistent quadrupole pairing interaction that recovers the Galilei invariance. 
The selfconsistent treatment of the pairing is important especially in the case of the density-dependent pairing 
force and in nuclei near drip- lines since in such cases the pairing potential A(r) in the surface region becomes 
relatively large. 

We have also checked that the selfconsistency in the pairing is fulfilled in the present calculation, by 
investigating the Nambu-Goldstone mode in the response functions for the monopole pair operators Pq = 
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/ drPt(r) and P = (Pj)*. When the HFB static pairing potential A(r) is not zero, the pairing rotational 
mode that has zero excitation energy should emerge as a Nambu-Goldstone mode associated with the nucleon 
number conservation of the Hamiltonian, or in other word the invariance with respect to the gauge rotation 
e l8N Q (N — J dr J2 a i/j^ {ra)ijj(ra) being the neutron or proton number operator ). The monopole response 
calculated for 24 exhibits the pairing rotation mode at the excitation energy very close to zero. It is found 
that we can move the energy of the pairing rotational mode exactly at zero energy just by modifying tinily 
the pairing force strength Vq within 1%. We also checked that the calculated monopole strength function 
for the nucleon number operators N carries no spurious strengths. This is again achieved by including 
selfconsistently the dynamical pairing correlation in the linear response equation. 



V. CONCLUSIONS 



We have formulated a new continuum linear response theory on the basis of the Hartree-Fock-Bogoliubov 
formalism in the coordinate space representation. This enables us to describe the pairing correlation in 
nuclei near drip-lines in a selfconsistent way both in the static ground state and in the dynamical collective 
responses. The formalism is able to include effects of the one-particle and the two-particle continuum states 
on the collective excitations. 

We have described the quadrupole response in the drip-line nucleus O with use of the density-dependent 
zero-range interactions. The strength functions for the quadrupole transition moments are obtained up to 
the giant resonance region. A low-lying mode which has significant neutron collectivity and is embedded 
in the neutron continuum states is obtained. We have analyzed in detail pairing effects important for the 
low-lying neutron mode. It is found that the low-lying neutron mode is sensitive to the density-dependence 
of the pairing correlation especially near the drip-line. The collective excitations in paired systems induce 
fluctuations not only in the normal density but also in the pair densities. The residual part of the pairing 
interaction causes dynamical correlation effects on the responses through the pair density fluctuations. The 
present theory describes selfconsistently both the static pairing effect caused by the HFB pair potential and 
the dynamical pairing correlation in the linear responses. The energy weighted sum rule is satisfied very 
accurately. This is because the selfconsistent treatment of the pairing restores the Galilei invariance, which 
would be violated if the HFB pair potential alone was included. 
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